A generalized Chudley-Elliott vibration-jump model 
in activated atom surface diffusion 
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Here the authors provide a generalized Chudley-Elliott expression for activated atom surface diffu- 
sion which takes into account the coupling between both low-frequency vibrational motion (namely, 
the frustrated translational modes) and diffusion. This expression is derived within the Gaussian ap- 
proximation framework for the intermediate scattering function at low coverage. Moreover, inelastic 
contributions (arising from creation and annihilation processes) to the full width at half maximum 
of the quasi-elastic peak are also obtained. 

PACS numbers: 



I. INTRODUCTION 

In 1960 Chudley and Elliott proposed a jump diffu- 
sion model for neutron scattering from a liquid, where 
the latter was assumed to be locally a quasi-crystalline 
structure. Within this model, the diffusive motion un- 
dergone by the (liquid) atoms is described in terms of 
large discrete jumps; otherwise, the atoms remain oscil- 
lating around their fixed (lattice) point between jumps. 
Because of the very different timescales ruling each type 
of motion (diffusive and oscillatory) , they are considered 
as completely uncorrelated. This jump model has also 
been extensively used to interpret the full width at half 
maximum (FWHM) of the so-called quasi-elastic peak 
{Q peak) observed in quasi-elastic helium atom (surface) 
scattering (QHAS) experiments ,Si2iii5i£ This peak, ob- 
tained from time-of-flight measurements converted to en- 
ergy transfers, accounts for the adatom surface diffusion 
process. Apart from the Q peak, in this type of exper- 
iments one also observes oscillatory motions associated 
to the temporary trapping of the adatom inside surface 
potential wells. This is the so-called frustrated transla- 
tional mode (or T mode), the lowest frequency mode of 
the adsorbate, whose lineshape is located very close to 
the Q peak which is centered at zero energy transfer. 

Following the Chudley-Elliot model, the coupling be- 
tween the diffusive and vibrational motions is very often 
neglected in theoretical descriptions. Of course, for short 
parallel momentum transfers of the scattered particles 
(He atoms), the coupling between both motions can be 
naturally disregarded, since no overlapping between the 
corresponding peaks will be observed. However, the same 
does not hold for large momentum transfers. This is im- 
portant because the experimental deconvolution of the Q 
and T-mode peaks is often carried out independently. 2,3 
In this paper, we propose a working formula for the whole 
lineshape based on the Gaussian approximation and that 
takes into account a generalized Chudley-Elliott expres- 
sion dealing with the coupling between the diffusive peak 
and the low-frequency vibrational modes. In this way, 
the experimental deconvolution procedure can be carried 



out in a more appropriate manner, and the typical infor- 
mation extracted (jump distributions, and diffusion and 
friction coefficients) from both motions at the same time 
would lead to more realistic values. 



II. THEORY 

In real diffusion experiments performed by means of 
the QHAS technique, the observable magnitude is the 
differential reflection probability, 7 



d 2 K(AK,Lu) 
dfldw 



= n d FS(AK,Lu) 



(1) 



This magnitude gives the probability that the probe (He) 
atoms scattered from the diffusing collective (chattered 
on the surface) reach a certain solid angle f2 with an en- 
ergy exchange fttu — Ef — Ei and parallel momentum 
transfer AK = K f - K 4 . In Eq. (TT|), n d is the (dif- 
fusing) particle concentration, and F is the atomic form 
factor, which depends on the interaction potential be- 
tween the probe atoms in the beam and the adparticles. 
The function S(AK, ui) is the so-called dynamic structure 
factor, which provides a complete information about the 
dynamics and structure of the adsorbed particle ensem- 
ble, and therefore also about the surface diffusion pro- 
cess and vibrational spectroscopy of their low frequency 
modes. Information about long distance correlations is 
obtained when using small values of AK, while long-time 
scale correlation information is available at small energy 
transfers, Hu>. 

In a series of previous works^££ we have analyzed 
in detail the lineshape structure of both the Q and T- 
mode peaks in terms of the standard jump Chudley-Elliot 
model. This model can be improved in order to provide 
a better analysis of the experimental lineshapes. For this 
goal, we start expressing the dynamic structure factor 
that appears in Eq. fl} as 



S(AK,c) = i J e-^I(AK,t) 



dt, 



(2) 
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where 

/(AK t) = ( e - iAK -[ R ( t )- R (°)]} = l e -i&KS*v Avi {t')dt'\ 

(3) 

is the so-called intermediate scattering function (the 
brackets denote ensemble averaging), with R(£) the tra- 
jectory of the adparticles and i>ak the corresponding ve- 
locity projected onto the direction of the parallel mo- 
mentum transfer AK (with length A A' = ||AK||). Both 
the dynamic structure factor and the intermediate scat- 
tering function can be readily obtained from standard 
Langevin numerical simulations.— Note that this type of 
simulation is valid as long as no collective, clustering, or 
self-organizing behaviors are present; in such cases one 
should make use of more detailed numerical techniques 
(e.g., molecular dynamics), which allow us to follow the 
processes that happen during short timescales (as the 
aforementioned ones). 

After a second order cumulant expansion in AK in the 
r.h.s. of the second equality of Eq. ([3]), I(AK, t) can be 
expressed as 

I(AK,t) « e -AK 2 /o'(t-*')C(«')dt' ; ( 4 ) 

where C(t) = (vak(0)uak(£)) is the velocity auto- 
correlation function. This is the so-called Gaussian 
approximation ^ 10 ^ 11 which is exact provided veloc- 
ity correlations at more than two different times are 
negligible: 9 - If non-Gaussian corrections are considered 
the theory cannot be expressed in a simple analytical 
manner; deviations from the Gaussian approximation are 
a clear indicator of the importance of such corrections. 

At low coverage, adsorbate-adsorbate interactions can 
be neglected; diffusion is then fully characterized by only 
studying the dynamics of an isolated adsorbate on the 
surface. This is the so-called single adsorbate approxi- 
mation in self-diffusion. Within the Langevin framework, 
the adsorbate-substrate interaction is described in terms 
of two contributions: (i) a deterministic, fitted adiabatic 
potential, V, which models the adsorbate-surface inter- 
action at zero surface temperature T s = K, and leading 
therefore to a deterministic force, F = — W; and (ii) 
a stochastic force (Gaussian white noise) accounting for 
the vibrational effects induced by the temperature on the 
(surface) lattice atoms that act on the adsorbate. 

For an almost flat surface (V fss 0), any direction is 
equivalent and therefore the dimensionality of the nu- 
merical Langevin simulation reduces to one. The corre- 
sponding numerical velocity autocorrelation function (in 
self-diffusion) then follows an exponential decreasing be- 
havior, 

C (t) = M± e -7t ) (5) 

m 

where m is the adparticle mass, ks the Boltzmann con- 
stant, and 7 the friction coefficient, which is related to 
the Gaussian white noise source through the fluctuation- 
dissipation theorem. Introducing ([5]) into Eq. Q, one 



obtains 10 

7(AA,t) = exp[- X 2 (e^ t + 7 i-l)], (6) 

where x is the shape parameter, defined asi 2 - 

X = VW) AK/j = lAK, (7) 

Here, W (v 2 ) — y^ksTs/m is the thermal velocity in one 
dimension, and I = y/ (v 2 ) /j is the mean free path. 
The diffusion coefficient is related to the friction as 
D = (v 2 )/~f. The only information about the structure 
of the lattice is found in the shape parameter [Eq. 10] 
through AK . When large parallel momentum transfers 
are under consideration, both the periodicity and struc- 
ture of the surface have to be taken into account. Con- 
sequently, the shape parameter should be changed for 
different lattices. 

The simplest model including the periodicity of the 
surface is that developed by Chudley and Elliott,— who 
proposed a master equation for the pair-distribution 
function in space and time^ assuming instantaneous dis- 
crete jumps on a two-dimensional Bravais lattice. The 
Fourier transform of this master equation gives rise to 
an exponential function for the intermediate scattering 
function, 

7(AK,t) = /(AK,0)e^( AK ' i . (8) 

The inverse of the correlation time, r„(AK), has a peri- 
odic dependence on AK of the form 

r„(AK) =^Pj[l-cos(j.AK)], (9) 
j 

where v is the total jump rate out of an adsorption site, 
and Pj is the relative probability that a jump with a 
displacement vector j occurs. The time Fourier trans- 
form of Eq. {5J is a Lorcnztian function with FWHM 
equal to r„(AK). At sufficiently small values of the par- 
allel momentum transfer, r„(AK) displays a cuadratic 
dependence on \ according to Eq. ([7|), 

r„(AK) « 2DAK 2 = 2 7X 2 . (10) 

From Eqs. ((9]) and (|10|) the shape parameter can be now 
generalized in terms of the jump distribution on the sur- 
face Bravais lattice, defining 

„(AK, S ^p) <n, 

which is valid for any value of the parallel momentum 
transfer, and where the subscript I denotes the impor- 
tance of the lattice. Note that, as a consequence of the 
effects induced by the periodicity and structure of the 
surface Bravais lattice, the shape parameter \ f° r f ree 
diffusion is replaced by Xl- This is related to the expo- 
nential dependence on time displayed by the intermediate 



3 



scattering function. Equation (fTTj) can then be seen as a 
definition for the new shape parameter, Xl, which takes 
into account the surface structure. 

Numerical Langevin simulations for corrugated surface 
potentials have shown^ that the Gaussian approximation 
is quite often appropriate. In such cases, the velocity au- 
tocorrelation function can be fitted to a general function 
described by 



C(t) = ("ak) e " 7 * cos(w T t + 5), 



(12) 



where lot gives the T-mode frequency (including the shift 
due to the surface temperature) and d is a phase which is 
added to better fit numerical results. The lot frequency is 
related to the harmonic one, loq, according to the relation 



LOT 



L0l 



11 

4 



(13) 



Analogously, in the case of a harmonic oscillator, 5 can 
also be expressed^ in terms of a relationship between 7 
and lot [in particular, 5 = (tan) _1 (7/2o;T)]- However, 
this constraint is lost when we apply Eq. (fl"2)) to cases of 
interest due to anharmonicities; as seen below, in such 
cases, 7, lot, and S are considered as fitting parameters. 

Introducing (TT2"]) into Eq. (TJJ the intermediate scatter- 
ing function can be expressed in terms of a double infinite 
sum of exponential functions, each contribution coming 
from the different quantum states of the T-mode oscilla- 
tor. This allows to use xi instead of x, an d, in analogy to 
the result obtained in Ref. 0, the intermediate scattering 
function of the full process can be generalized as 

7(K,<) = cxp [-xU(ut,5M 

fA 1 - x fA 2 t (~l)"(-l) m 2(n+m) 



e -^ Al - x - A2t v "1,1,"' XT^'AZA 

n.m 



n\m\ w " 3 4 



where 



Ai 



An 



7 2 [2jlot sin S + (lo t - j 2 ) cos S] 
(7 2 +^) 2 : 
7 2 (7 cos 6 — lot sin S) 



■y 2 +LO T 



(14) 

(15a) 
(15b) 
(15c) 

„ . vr ( 15d ) 

2(7 + IL0 T ) 

Now, Fourier transforming Eq. (|14[) . an analytical ex- 
pression for the dynamic structure factor is also readily 
obtained, 



A, = 



Aa 



7 



2(7 — ilot) 2 ' 

,2 



7 



S(AK,w) 

e 



TT 



-i(m— n)5 



r I m I 



n.m— 



[ X 2 A 2 + (n + mfr] 



[cj - (n - m)w T ] 2 + [% 2 A 2 + (n + 771)7] 5 



(16) 



where the n and m indexes represent different excitations 
of T-mode creation and annihilation events. 

The double sum over Lorentzian shapes in Eq. (fl6|) 
clearly shows that both motions (diffusion and oscilla- 
tions) cannot be separated at all. Only the terms with 
m = n (events with a zero energy balance) contribute 
to the Q peak, while those with n 7^ m will contribute 
to the T-mode peaks. As discussed in Refs. and [H, 
the shape parameter also governs the shape of the peaks. 
Depending on the new shape parameter (the dependence 
on AK is only through xi), the sums will contribute 
globally to different intermediate shapes ranging from a 
Gaussian function (high values of xi) to a Lorentzian one 
(Xi <§C 1). This smooth shape transition is the so-called 
motional narrowing effect, well-known in the theory of 
nuclear magnetic resonance lineshapes^ For a given sys- 
tem, this effect could be observable by changing the AK 
value, as shown below. In Ref. we proposed a similar 
functional form to Eq. (|16[) as a working formula, but 
using instead the shape parameter for a flat surface [x, 
as given by Eq. ([7])]. With the new expression provided 
here, where x has been replaced by Xh one can now de- 
termine the lineshapes displayed by the Q and T-mode 
peaks in numerical calculations or in experiments after 
deconvolution for the whole first Brillouin zone. This 
is an important issue because the experimental deconvo- 
lution procedure is generally carried out for each peak 
separately (i.e., there is no overlapping), pre-assuming a 
given shape (an effective Lorenztian function) with no 
theoretical justification. This assumption leads to a pro- 
portionality between the effective width and the diffusion 
coefficient, which is not correct in general. In this sense, 
to perform the whole deconvolution procedure a gener- 
alization of such an equation is convenient, which will 
allow to extract, at the same time, good diffusion coeffi- 
cients, information about jump mechanisms, friction co- 
efficients, and low T-mode frequencies. When more than 
one Lorentzian function contribute to Eq. (| 16[) , diffusion 
coefficients should be then obtained from the calculated 
jump distributions, and therefore not directly extracted 
from the experimental Q peak. 

From the generalized expression for the lineshape of 
the full process, given by Eq. (fTB"]) . the width of each sum 
contributing to the Q peak reads now as 



r£(AK) = 



r„(AK) 

27 



A 2 + 2nj. 



(17) 



Note that A2 provides information about the dynami- 
cal process through the frequency associated to the T- 
mode, the friction coefficient and the phase of the ve- 
locity autocorrelation function [see Eqs. (fl"2"]) and (|15bj) ]. 
The surface temperature dependence appears implicitly 
in the shift of the position of the T-mode frequency and 
the jump distribution Pj through r„(AK) according to 
Eq. ([5]). Typically, the most important contribution to 
the Q peak at low temperatures is r^ =0 . Within this 
model, if 5 = 0, lot = 0, and the regime of very small val- 
ues of the parallel momentum transfer is studied, Eq. |TT 



4 



is then satisfied and T^(AK) = DAK 2 , as obtained 
when a pure (or effective) Lorentzian function is assumed. 
Inelastic contributions (with zero energy balance), r^ , 
to the total width due to creation and annihilation pro- 
cesses of the T-mode are important with the surface tem- 
perature and increases by integer numbers of 7. Mixed 
Lorentzian-Gaussian lineshapes are very often obtained, 
the FWHM being calculated from the general expression 
of the dynamic structure factor. 

In a similar fashion, the width of each sum contributing 
to the T-peaks, according to Eq. IT6|) . is given by 



rI +m (AK) = 



r„(AK) 

27 



A 2 + {n + m) 1 , (18) 



where the jump distribution Pj is also involved in the vi- 
brational relaxation process through r t/ (AK). It is gen- 
erally assumed that the corresponding width is 7, which 
(within this model) is obtained when n = 1, m = 0, 
6 = 0, and 7 <C <^o; note that the first term in Eq. (fT8|) 
behaves as j 2 /ujq, which is much smaller than the second 
term (which is equal to 7). In this case, the lineshape is 
a pure Lorentzian function. 



III. RESULTS 

As an illustration of the theoretical approach proposed 
and described in Sec. |TTJ we have carried out numer- 
ical Langevin simulations for the Na/Cu(001) system, 
which has been largely analyzed both experimentally 
and theoretically The adiabatic nonscparable adsorbate- 
substrate interaction is taken from Ref. U with a cover- 
age of 2.8% and for a surface temperature T s = 150 K. 
In Fig. [IJa), the shape parameter \i is plotted as a 
function of the parallel momentum transfer covering the 
first Brillouin zone. The solid curve has been obtained 
from Eqs. ((9]) and (jTTJ) after jump distributions were 
calculated. The circles give the fittings of the differ- 
ent Langevin simulations to the dynamic structure factor 
given by Eq. (fTo)) . In panels (b) and (c) of Fig. [TJ the 
quality of the fitting for AK =1.23 A -1 [black circle in 
panel (a)], for the Q and the two T- modes peaks, respec- 
tively, are shown; dotted lines correspond to the applica- 
tion of Eq. (|16p with no adjustable parameters (nominal 
values: 7 = 0.6 meV, lot — 5.99 meV, and 6 = 0.05), and 
dashed lines to the case where Eq. (fTB) has been used as 
a working formula with three different adjustable param- 
eters (7 = 0.53 meV, uj t = 5.82 meV, and 6 = 0.032). 
The whole spectrum is better reproduced in the last case; 
in the first fitting the right positions of the T-mode peaks 
are obtained but with very small intensities compared to 
the numerical simulations. Note, however, that for the 
Q peak, both fittings are approximately the same. The 
small deviations from the nominal values of the three 
parameters are a clear indicator of the non-Gaussian be- 
havior of the dynamical system. Concerning the shape 
of the Q peak, when \i is greater than one (central re- 
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FIG. 1: (Color online.) (a) Shape parameter \i as a function 
of the parallel momentum transfer covering the first Brillouin 
zone for the Na/Cu(001) system at T 3 = 150 K. Open circles 
correspond to the fitting of the numerical simulation, while 
the solid line is the result obtained from Eq. (|11[) . The dy- 
namic structure factor for the maximum momentum transfer, 
AK = 1.23 A" 1 (solid circle), is plotted in panels for the Q 
peak (b) and for the two T-mode peaks (c). Solid (black) lines 
corresponds to the numerical simulation; dotted (blue) lines 
represent the result obtained from using the working formula 
with no fitting parameters; dashed (red) lines are the result 
when three fitting parameters (7, ujt, and S) are employed. 
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gion of the Brillouin zone) , a mixed Gaussian-Lorentzian 
function is expected. 

In Fig. [21 the convolution of the dynamic structure 
factor (denoted by Cs) obtained from the Langevin nu- 
merical simulation for AK = 1.23 A" 1 [see Fig. [He)] i s 
plotted (solid line). Moreover, the convolution of sepa- 
rate Lorentzian functions for the Q and T-mode peaks 
are also shown (dashed line for the Q peak and dotted 
lines for the T-mode peaks). The two panels correspond 
to two different surface temperatures: (a) 150 K and (b) 
250 K. The response function of the experimental setup 
is known and given in Ref. 0. Discrepancies between the 
standard procedure and the one proposed here are appar- 
ent, as indicated by their different backgrounds, widths, 
and peak intensities. 

As a word of conclusion, we consider that a convolution 
of the whole spectrum should be carried out in order to 
properly take into account the contributions of the events 
with zero energy balance in the Q peak, as well as the 
existing overlapping between both Q and T-mode peaks 
at larger momentum transfers and high surface temper- 
atures. 
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FIG. 2: (Color online.) Numerical convolution (solid line) of 
the Q and T-mode peaks for the Na/Cu(001) system with 
parallel momentum transfer AK = 1.23 A -1 at two different 
surface temperatures: (a) 150 K and (b) 250 K. Dashed (red) 
lines correspond to the convolution of a pure Lorentzian func- 
tion for the Q peak; dotted (blue) lines denote the same, but 
for the T-mode peaks. 
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